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ABSTRACT 


This thesis addresses the problem of integrated design of the aircraft plant pa- 
rameters and of the corresponding feedback controller. The plant parameters are 
typically the sizes of the control surfaces or other aerodynamical surfaces of the 
aircraft. The approach is to rewrite the aircraft dynamic requirements as linear 
matrix inequalities (LMI’s) and to optimize a linear cost faction associated with 
aircraft plant parameters, while meeting the LMI constraints. An algorithm using 
Matlab and LMI-Lab has been developed. This algorithm has been used for inte- 
grated plant/controller design of an F-4 aircraft at five different flight conditions. 
The result consists of a set of controllers, one for each flight condition, and a single 


solution for the optimal sizes of the F-4 stabilator and spoiler control surfaces. 
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I. INTRODUCTION 


Control systems of modern aircraft have become increasingly complex. That 
is especially true in the case of fighter aircraft with relaxed stability and multiple 
control surfaces working in each axis. 

The control system, control surfaces included, must satisfy all the performance 
requirements (stability, disturbance rejection, maneuverability etc.), while not causing 
excessive drag, weight or signatures. It is thus interesting to study the problem of 
finding the smallest size of some parameters describing the physical configuration of 
the aircraft, while still being able to find a controller that satisfies the performance 
requirements. 

While the importance of integrated plant and controller design has been recog- 
nized, little had been done in the area prior to [Ref. 1]. In this thesis the ideas of 
integrating a minimization of the control surfaces with the design of a controller, that 
were developed in [Ref. 1], are developed further. 

The design procedure is based on a set of requirements that are to be fulfilled at 
a number of different flight conditions and a cost function for the plant parameters, 
e.g. the size of the control surfaces. The result consists of a set of controllers, one 
for each flight condition, and one single optimal solution for the size of the plant 
parameters. The eouoaad design procedure is still in an early prototype phase, but 
shows promise of becoming a useful tool. 

The thesis is organized as follows. Chapter II contains the formulation of the 
treated problem as an H,, control problem, and the parametrization of the plant. In 
Chapter III the algorithm from Chapter IT is applied to the longitudinal dynamics 


of an F-4 fighter aircraft. The plant parameters are chosen to be the sizes of the 
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moving horizontal tail and of the spoiler. Chapter II also contains the results from 
the plant optimization and the associated controllers. The program code in Matlab 


and LMI-Lab that is used to solve this problem is included in the appendix. 











Il. PROBLEM STATEMENT 


The task of designing the control system of an aircraft has become increasingly 
complex in recent time. That is especially true in the case of modern fighter aircraft 
with relaxed stability and many control surfaces working in each axis. The control 
system must, of course, meet all the imposed performance requirements, but as ex- 
cessive control capability costs weight, space, money, signature etc., it is desireable 
to optimize the control system. 

The problem to be addressed is to find the minimal size of some plant parameters 
(typically the size of the control surfaces or other aerodynamical surfaces), while 
finding a feedback controller that meets the stated control requirements. The design 
requirements are given for several flight conditions and the designed combination of 
plant parameter and controller must, of course, work for all of the flight conditions. 
The result of the optimization is one set of values for the plant parameters and a 
set of controllers, one for each flight condition. The objective is to use this set of 
controllers to design a gain scheduled controller, that works for all the given flight 
conditions (and hopefully in between those given conditions). 

When the sizes of the surfaces are varied, the dynamics of the aircraft changes. 
It is therefore natural to integrate the minimization of a cost function associated with 
the plant parameters with the controller design. The controller is designed using Hq 
methods, which are summarized in Section A. The problem of minimizing the plant 


parameters, while still being able to find a feasible controller, is treated in Section B. 














A. H,, CONTROL 





Figure 2.1: Feedback interconnection of plant G and controller K. 


A general feedback system with exogenous input w and output to be controlled 
z is shown in Figure 2.1. In this context we will only consider linear, time invariant 
systems. Here w is a reference signal that includes disturbance inputs. We wish 
the system to reject the disturbance and track the reference signal, i.e. we want 
a controller that makes the difference between the reference input and the actual 
output as small as possible, while not expending too much control energy. This 
is accomplished by making z, which is a combination of tracking errors and the 
control inputs, as small as possible. Let T,.,(G,K) denote the closed loop transfer 
fuction from the input w to the output z. The solution to the H,. control problem 
is the controller K that minimizes the infinity norm of the transfer function, i.e. 
ToC, K) 


The infinity norm is a generalization of the maximal gain and is defined as the 








took! 


supremum over all frequencies of the maximal singular value of the transfer function, 


i.e. 


|Toollso = sup{F(Tow(ju))} 


Although the above definition is given in terms of the frequency domain, state- 


space methods are very useful in H,. synthesis. A state-space realization of the plant 











depicted in Figure 2.1 with state feedback can be written as 


z = Cx+Du (2.2) 


Where xz € R",z € RP, ue R™ and we R!. 

To ensure the existence of a solution to the H,, control problem, we have to 
assume that (C2, A, Bg) is stabilizable and detectable, and dee D has fitiéaaly inde- 
pendant columns. The inputs and outputs, i.e. w and z, are usually scaled in such a 
way that the design requirements are satisfied if ||Tzu(G, K)|lo < 1. 

It can be shown [Ref. 2] that there exists a controller K such that ||T(G, K)|loo < 
1 if and only if there exist W and Y, with Y symmetric and positive definite, such 
that 


AY+YAT+BW+4+W7B? B, (CY+DW)? 
R(W,Y) = Br —I 0 <0 (2.4) 
CY + DW 0 —I 


If the linear (or rather affine) matrix inequality (2.4) has a solution, the con- 


troller K = WY! solves the H,. control problem. 


B. PLANT PARAMETRIZATION 

We consider a dynamical system at several given operating conditions, or more 
specifically an aircraft with performance requirements given for several flight condi- 
tions. At each of the flight conditions we use a linear model to describe the dynamics 
of the aircraft. The linear models are obtained by linearizing the general non-linear 
six degree of freedom equations of motion for each flight condition. The linearized 
models are different for each flight condition and with 7, 7 = 1,2,...1, denoting the 


respective flight conditions we have 














¢ = A'e+ Biw+ Biu (2.5) 


z = C'r+ Du (2.6) 


Where z € R® is the state-vector representing perturbations around the trimmed 
condition, w € R% is the external input, u € R™ is the control input and z € RP is 
the ouput to be regulated. 

We now parametrize the plant with the plant parameter ¢ = iis CorcsenCeles 
where ¢ is normalized such that 0 < ¢; < 1. Associated with ¢ is.a linear cost function 
= c?(, which we wish to minimize. The A’, Bi and Bi, matrices, describing the 


plant are parametrized as the first terms of the Taylor expansion around the point 


¢=0as 
k 

AY(¢) = Ab + 2 AG (2.7) 

Bi(¢) = Bla + BG (2.8) 

BQ) = Bis + BG (2.9 

where 

A= Fl. (2.10) 

B= Fel. (2.11) 

B3; = Sel... (2.12) 


The Ci and D' matrices could also depend on ¢, but that case is not considered 


here. 


For each of the flight conditions we can find a steady state H.. controller by 


solving (2.4) for W* and Y*. However, the plant matrices now depend on ¢, hence 








i ANC)Y? + Y'A'T(C) + B(C)W! + WT BIT(C) BY(G) (C*Y! + Diw')t 
R(W',Y',¢) = Bi" (¢) =f 0 <0 
c'y'+ Dwi 0 ~I 


Our problem is thus 
min J = c¢ 
s.t. (Y',W',C) €® 
where the feasible set ® is defined as 


= {Wie R™ YER C:Y'>0,Yi=V'70<6G <1, R(WY',C) <0} 


C. ALGORITHM 

Following the results in [Ref. 1] the optimization problem in W’, Y* and ¢ is 
solved as follows. We observe that R, for a fixed ¢, is affine in the two remaining 
variables W‘ and Y'. We further observe that R, for fixed Y' and W’, is affine in the 
sole remaining variable ¢. 


Hence the approach is to: 
1. Fix ¢ = ( and find W' and Y* > 0,2 = 1,2,...,/ such that R(W',Y', 00) <0 


2. With W* and Y' from step 1, find the ¢ with ¢; > 0 that minimizes J = cT™¢ 
s.t. R(W',Y',C) <0 


3. Go to step 1. until exit criterion is satisfied. 


The problem to be solved in step 1 is to solve 1 LMI’s for the matrices W? 
and Y'. The / LMI’s are, however, independent of each other and can be solved in 
sequence. The problem in step 2 is that of solving one single LMI for the minimal 
plant parameter ¢. Here all flight conditions must be treated simultaneously, as the 


same physical aircraft must be used at all the flight conditions. 


vs 














The solution of LMI’s is supported by software packages such as LMI-Lab [Ref. 
6], which is used in the example given in Chapter II]. While the approach to formu- 
late the problem as two LMI’s solved in sequence makes the calculations manageble, 
there is unfortunately no guarantee that the resulting algorithm will converge to a 
stationary point, [Ref. 3]. Neither is it clear whether or not the feasible set © is 
convex. Convexity is, of course, of utmost importance in every optimization problem, 
as when a convex function is minimized over a convex set, all minimas are global. 

Chapter III describes the application of this algorithm to the longitudinal dy- 


namics of an F-4 fighter aircraft. 








Ill. LONGITUDINAL CONTROL OF F-4 


The specific example concerns the control of the longitudinal dynamics of an 
F-4C fighter aircraft. The objective is to size the two control surfaces, elevator and 
spoiler, and to design a controller for elevator, spoiler and thrust for each of the 
different flight conditions. Five flight conditions, three subsonic and two supersonic, 
with data from [Ref. 4] are used; see tables 3.1 and 3.2. . . 

The problem formulation is basically the same as in (Ref. 1], i.e. to minimize 
the sizes of the elevator and the spoiler, while finding a controller that satisfies a 
set of performance requirements for each flight condition. The main difference is 
that thrust control is included and that multiple flight conditions are treated. The 
numerical algorithm, which uses LMI-Lab, is completely different from the one used 
in [Ref. 1]. 

In section A the equations of motion are parametrized according to the chosen 
parameters, elevator and spoiler. Section B treats the performance requirements and 
how they are expressed in terms of H,, control. Sections C-E contain the results from 
the performed calculations. 


Nomenclature: 


X, Aerodynamic derivative with respect to states or controls in X-direction 
Z, Aerodynamic derivative with respect to states or controls in Z-direction 
M, Aerodynamic derivative with respect to states or controls in M-direction 
I, Moment of inertia around the y-axis 

m Mass 

U Trimmed velocity in x-direction 














ao Trimmed angle of attack 

1 Distance from center of gravity to elevator 
6. Elevator deflection 

6, Spoiler deflection 

6; Thrust perturbation from equilibrium 

0, Disturbance in angle of attack, a 


o., Disturbance in flight path angle, 7 


A. EQUATIONS OF MOTION AND AERODYNAMICAL DATA 

The equations describing the dynamics of an aircraft in six degrees of freedom 
consist of a set of coupled non-linear differential equations. These equations are 
conveniently formulated in a body fixed coordinate system. In the design of the 
control system, we generally consider the set of linearized equations describing the 
dynamics around a certain nominal flight condition. 

In this section we use a model with four linear differential equations to describe 
the longitudinal dynamics. In general the longitudinal dynamics of an aircraft can be 
divided into two modes: the short period and the long period. From [Ref. 4] we obtain 
aerodynamical data for an F-4 fighter aircraft operating at five flight conditions: three 
subsonic and two supersonic. 

The plant parameters to be optimized are the size of the moving horizontal 
tail and the size of a spoiler. These sizes are represented by the two coefficients 
Zs and Z,, which give the vertical acceleration per radian of elevator and spoiler 
deflection respectively. Hence, the plant parameter is ¢ = Caer peer ie The 
aerodynamical data is adapted to fit that parametrization, by first calculating the 
aerodynamical derivatives for the condition ¢ = 0 and then calculating all the matrices 


needed in the equations of motion. 
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TABLE 3.1: Aircraft Data (Ref. 4) 


[Fight Condition OM. ] eo] @] oO 










3 
[Altitude (f#) | 15,000-| 35,000 | 35,000 | 35,000 | 45,000_ 
[Mach | 09 [ 06 | 09 | 12 | 15 | 
aire) 189 
[-~m(slugsy | 1210 | i910 [ 1210 ra10_ | 7310_| 
[_T(slug FP) _| 24970 | 27360_[ 25040 [24970 | 25040 
3033 “3033 
26 

















The spoiler is an add-on to the data from [Ref. 4] and is modeled in the same 
manner as in [Ref. 1], i.e. as giving a maximal acceleration of 0.2g in the z-direction 
and a corresponding acceleration in the z-direction of —0.02g. 

We need to find out how the varying of the plant parameters affect the aerody- 
namical derivatives describing the aircraft dynamics. Changing the size of the moving 
tail affects the dynamics of the aircraft drastically. As the tail is reduced the stability 
as well as the damping is reduced. The other plant parameter, the spoiler, does, 
however, not affect any coefficients other than X, and Z,. The division of the aero- 
dynamic derivatives into one part that depends on the wing and body of the aircrft 


and one part that depends on the tail is done according to [Ref. 5] as 


a (3.1) 
I 

Zeus = Dy 7s (3.2) 
i 

Lith, = Za — 74s (3.3) 
1 

Mii. = “MyeeeZ5 (3.4) 
I, 
[? 

Maus = M,~—=Zs (3.5) 














TABLE 3.2: Nominal Dimensional Derivatives (Ref. 4) 


KotH) 
i170 

Zt file 
Zaft 
Z(fi/s) | 6.00} 1.82 | 2.89 | 4.09] 2.24 J 
Xysy_| 0.00 | 001] 0.00 [001 | 0.00 
ATED 
M;(s~*) : -I1. ; 

2 

0.919 | 0. 





——s 























—M(s-7) [25.00 [4.90 [11.40 | 20.70" | 16.00 
ZC fils4) 








Maw = Mg — Zs 


(3.6) 


The resulting wing/body aerodynamical derivatives are shown in Table 3.3. A 


for concern. 


comparison with the values in Table 3.2 shows that the removal of the horizontal tail 
causes the static stability to decrease. In the three subsonic conditions the aircraft 
is now unstable (M, > 0). The damping has also decreased, i.e. My + M, is now 
closer to zero. The fact that M, is positive is unusual, but as it is usually difficult 


to separate M, and M,, and they still add up to a negative number, there is no need 


For the longitudinal dynamics we use the linearized equations of motion [Ref. 4] 
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with the four state variables given as [u/U,a,q,6]?. These states represent perturba- 


tions around the trimmed values [U, ao, 0, 40]7. By using the dimensional aerodynamic 








TABLE 3.3: Wing/Body Dimensional Derivatives 












[Fittadinl Oe] e >a] o] 
[Xue ") [W025 [0.0176 [-0.0193 | 0.0136 | -0.0072 | 
[xefs*) [6.93 _[ 2625 [743] 16.83 | 20.98 | 
2s) [-i70 | -0337| sr | 077 [0505 | 
[_2.(Ft]") | -996.2 | -14l.2 | 495.15} “757.9 | -646.0_| 
| _Z.(Ft/s) | 0.652] 0.256 | 0.318 | 0.551 | 0.504_| 
[4tf/s) | -335_[ 097s [156 | 230 | -113_| 
[M(F-'s)_| 0.0001 0.0022 | 0.0025 | 
Mi) | 800 | 207 | 3s2 | -833 | -12.01 | 
M(t) [0.163 0.057 | o07T | 0.122 [0.130_| 
Ms) [0.375 | -0.109-[ -0.182 [0.386 | -0.256 | 





derivatives the equations of motion for the longitudinal case assume the form 
Tit = Aix + Bilw+ Bi,u (37) 


where w = [Qref,Yres]? is the external input and u = [6,6,,6;]? is the commanded 
values of elevator deflection, spoiler deflection and thrust. 
With aerodynamical data for the respective flight condition, the plant matrices 


are given as 














osx 
oo 
oocr 


be 
u=| 6, 
bt 


We now write the A’, B‘, and Bi matrices as 
3 1? 2 


2 
ANC) = Agt SO AiG 


j=l 


Bi(¢) = Biot > BijG 


j=l 


Bi(¢) = Bry + D7 Baj6j 


j=l 


(3.11) 


(3.12) 


(3.13) 


(3.14) 


(3.15) 
(3.16) 


(3.17) 


This is the equivalent of taking the first two terms in the Taylor series expansion 


around the point ¢ = [0,0]7. Using the matrices defined in (3.8-3.11) we have 


Ay = (hi) AS 














; . OA, OAb " OL )-* 3 in OAL 
Ay = re ~ O25 x5 Zénom =i O%, At, + (Ii)-? OZs ——)Zsnom 
» _ OA, adi = fOE Den eetOAG 
Ay = Ob = aZ, oar Lanpan = (—>- OZ, =a A, ee aZ, =7- )Zsnom 
Big o (Ti) Bin 
i _ OBig _ OBio pO 3 i\-1 9B} 
Buy -s OG, OZ; Zénom oe ( Oz: Bi, + +(I,)7 OZ; 
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*)Zénom 


(3.18) 
(3.19) 
(3.20) 


(3.21) 


(3.22) 


re 


t 
By 


needed derivatives. 


O(In)- 
OLs 








OB}, - dBi 
C2 OZ, 
(In) Ban 

OBi, = OB; 




















A  OZs 
OBi, 2 OB}, 
OC. 7 OZ, 


1 


snom — ( 
énom = ( 
snom — ( 





On) 
OZ, 





OI)” 
OZ; 
O(In)™ 


OZ, 








-1 
Bin + ( 


Bi, + ( 


B3,, + ( 


T 


nr 


r 


nr 


rT’ 


n 


) 


) 
) 





-1 OBin 


pz, Zenon 


OB; 
OZs 
-1 OB,, 
OZ, 


-1 


Veron 





Zenon 


(3.23) 
(3.24) 


(3.25) 


(3.26) 


Where all the matrices are evaluated at [Z;, Z,] = [0,0]. We now calculate the 








OZ, 





| 0 
aBi, | 0 
0 


Coo & 


(3.34) 


0 0 
B. DESIGN OBJECTIVES AND AUGMENTED SYSTEM 


The design objectives follow those outlined in [Ref. 1]: 


e Step response - There should be no steady state error in the response to a step 


input in angle of attack, a, or flight path angle, y = 6 — a. 
e Closed loop must be internally stable. 


e Gust disturbaces with a magnitude of 5 ft/s in the z-direction should not yield 
elevator deflection greater than 20 degrees, spoiler deflection greater than 40 


degrees or errors in angle of attack greater than 1.5 degrees. 


e A flight path angle of 1 degree must be generated using elevator and spoiler 


deflection of less than 20 and 40 degrees, respectively. 


To accommodate the steady state tracking requirement, the plant is augmented 
with the integrals of a and 7. The augmented system has six states, x = [u/U,a,q, 9, a/s,7/s]*. 


The parametrized augmented system is thus 


a 


2 2 2 
(Abo + > Au jGi)@a + (Baro + Ds Bir Sj)w + (Bizo + S By2;6;)u (3-35) 
j=l j=l jt 


z = Cir,+ Diu (3.36) 


with the respective matrices given as 


(3.37) 


ooooo & 
ooo oo © 


























Bix = (3.46) 
0 0 O 
0 0 O 
0 0 00 0 0 
0 0 00 0 0 
;_|0 0 00 0 0 
C= 0 a 0 0 0 0 (3.47) 
0 0 00 we 0 
00000 w 
— , 0 
0 Ssmaz 4 
Di = Bimaz 3.4 
0 0 0 
0 0 0 
be 
bs 
z=| 0 (3.49) 
a/s 
7/8 


The scaling of the in- and output to get ||T2w||o.. < 1 shows up through the og 


and o., in the By; matrices and through the weights in the C, and D, matrices. 


C. CALCULATIONS 

The calculations were done in Matlab. The code is listed in the appendix. The 
linear matrix inequalities are solved by using LMI-Lab, [Ref. 6], which is a recently 
developed toolbox to Matlab. W' and Y‘, i.e. the controllers, are obtained from the 
function ”feasp” and the minimization of ¢ is done by using ”linobj”. The function 


*feasp” sometimes, however, returns a solution that is not at all negative definite. 
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That problem is solved by testing the calculated solution and, when it is infeasible, 
discarding it in favor of the solution from the previous iteration. When ”feasp” is 
again called in the next iteration, it usually gives a feasible solution. 

The code is fairly fast - a run takes about 10 minutes on a Sparc 10. 


The calculations were performed with the following values: 


G,= 5/U; 
bemar = 20 degrees 
Osmax = 40 degrees 


Dine = 0.05g 


D. OPTIMIZATION RESULTS 

As we are uncertain about whether or not the feasible set © is convex, the con- 
vergence properties of the algorithm were tested to see if the solution would converge 
to a single optimal point or not. 

The convergence properties were tested by running the algorithm with different 
initial conditions. First we varied the initial value of the elevator, i.e. Cyinitial, between 
0.1 and 1.5, while holding the initial value of the spoiler as the nominal, i.e. Coinitiol = 
1. In 3.1 we see the optmiaation history for each iteration of the two plant parameters. 
The solution tends to converge to smaller values, the smaller the initial value is. The 
same results are shown in 3.2, where the two values ¢, and (2 are plotted against each 
other. The cost function for the different initial values of ¢; is shown in Figure 3.3. 
It is obvious that the algorithm does not converge to the same solution, when it is 


run with different initial values. 








The initial value of the spoiler, ie (2, was varied between 1.5 and 0.6, while 
the initial value of the elevator was held constant at 0.5. The optimization history 
is shown in Figure 3.4, while Figure 3.5 shows ¢, as function of ¢;. In this case the 
algorithm converges to the solution ¢,,: = [0.06,0.55]7, for all initial conditions but 
the two smallest. Figure 3.6 shows the cost function for different initial values of 2. 

The conclusion is that the calculated solution is sensitive to the initial value, 


especially the initial value of ¢,. The obvious conjecture is that the feasible set ® not 


is convex. 





Figure 3.1: Optimization history of ¢ - varying C1initia! 


The calculated solution always gives a very small size of the elevator. That is, 


however, consistent with the problem formulation. 


E. CONTROLLERS 
The controllers for the five flight conditions that are listed below, were calculated 
for ¢ = [0.1,0.7]" which is close to the lowest values given in Figures (3.1-3.5). The 


resulting five controllers are presented in below. Tables (3.4-3.8) show the eigenvalues 
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“o 05 1 1.5 


Figure 3.2: Optimization history: (2 as function of ¢ 


and dampings for the closed loop system for the five flight conditions. The designed 
controllers give the closed loop system rather slow dynamics. This can be corrected 
by adjusting the weights in the C matrix. The damping, which varies between 0.54 
and 0.88 is, however, good. As the purpose of the performed work is to develop a 
procedure for the design of an controller, rather than actually designing a controller, 


the results are clearly acceptable. 








39.7. 10.7 42 29.9 26.0 10.0 

K'=| 505.3 -149.3 0.9 150.5 —130.5 32.7 
—~1201.0 297.1 —1.7 —300.7 258.9 —65.1 

—~3.76 36.86 17.07 5.77 57.94 2.70 

kK? = 105.6 —148.4 -2.04 139.9 -—39.57 31.30 
—154.97 71.75 -—0.246 —72.83 4.70 —17.23 

15.28 -2.55 7.49 27.44 40.65 ~~ 7.74 

Ke =| 230.2 —151.4 0.913 151.9 —59.50 32.15 
—362.5 153.8 —1.101 -156.1 55.94 —31.67 


459.4 —213.9 1.06 215.3 —100.1 41.96 


38.51 —19.58 3.80 33.03 29.71 8.88 
Kkt= 
—837.7 295.6 —1.44 -—298.4 134.5 —57.38 
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Figure 3.3: Costfunction for different Qiinitig) as function of ¢ 


Ki =] 787.0 —505.4 2.4 510.5 -—149.9 89.2 


105.4 -—64.3 68 87.6 23.3 17.5 
~—1065.3 543.6 -2.6 -—550.0 153.4 —95.0 


TABLE 3.4: Closed loop eigenvalues and damping, condition 1 


Freq. (rad]sec) 















Figures (3.7-3.9) show the actuator step response for flight condition 1, while 


figures (3.10-3.12) show the Bode plots for the actuator dynamics. 
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0 5 10 15 20 25 30 35 


Figure 3.4: Optimization history of ¢ - varying Coinitial 


TABLE 3.5: Closed loop eigenvalues and damping, condition 2 


Freq. (rade) 










1.0000 





TABLE 3.6: Closed loop eigenvalues and damping, condition 3 


Freq, (rad/sec) 
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0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 


Figure 3.5: Optimization history: (2 as function of ¢ 





0 5 10 16 20 25 30 35 


Figure 3.6: Costfunction for different Cinitiai 
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TABLE 3.7: Closed loop eigenvalues and damping, condition 4 


[__Figenvalue [Damping | Freq. (rad[sec) | 
[-O1972 + 0.1doe | 0.8026 | O2aS | 
[0.1972 - 0.14661 | 0.8026 [0.2458 | 
[1.6795 | 1.0000 [1.6705 


[2.1119 [1.0000 | 2.1119] 
[37604 + 5.8035 | 0.5388 [6.9959 | 


[3.7694 - 5.80351 | 0.5388 | 6.9959 | 














TABLE 3.8: Closed loop eigenvalues and damping, condition 5 


[Eigenvalue | Damping | Freq. (rad/sec) J 
[607 + ore] 0.8330 | 0.2088 ‘| 
[0.1607 0.1127; [~0.8330_| 0.2038 
-_-20066 | 1.0000 | 9.0966] 


-3.4113 1.0000 [3.4113 | 
4/3355 + 3.93291 | 0.7407 5.8536 
-4,3355 - 3.93291 | 0.7407 5.8536 
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Figure 3.7: Actuator step response 
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Figure 3.8: Actuator step response 
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Figure 3.9: Actuator step response 
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Figure 3.10: Bode plot for elevator 





Input 2 


HKTOF reese ene EL ce dB dobeld Slee eeeede 
fae) 
uo 
© 20h 8 edocs eer 
oa 
S 

~30 st staves dea a tiars 

-40 : 

107 10° 10° 10' ° 
Frequency (rad/sec) 


10 





10” 10° 10 
Frequency (rad/sec) 


Figure 3.11: Bode plot for spoiler 
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Figure 3.12: 
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Bode plot for engine thrust 
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IV. CONCLUSIONS AND 
RECOMMENDATIONS 


A. CONCLUSIONS 

This thesis develops the technique of integrating the design of aircraft plant 
parameters and the design of the corresponding feedback controller. The developed 
algorithm using Matlab and LMI-Lab gives a solution to the problem of minimizing 
the control surfaces of an aircraft, while finding a feasible H,. controller. The result, 
when this algorithm is applied to the longitudinal dynamics of an F-4 aircraft, is a 
set of controllers, one for each flight condition, and a single value of the sizes of the 
moving horizontal tail and the spoiler. 

The approach used to formulate an integrated plant optimization and H,, con- 
troller design as a linear matrix inequality results in a non-convex optimization prob- 
lem. 

Compared to [Ref. 1] the problem formulation has been expanded to handle 
multiple flight conditions and to include thrust control. As the developed algorithm 
using LMI-Lab is more computationally efficient than the previous algorithm, it is 
possible to handle much larger problems than in [Ref. 1]. 

The calculated solutions do, however, depend on the initial value, as the problem 
is non-convex. It is thus nepeaaey to run the algorithm several times, to be reasonably 
sure that a good quality solution has been obtained. 

This work has increased the knowledge on how integrated plant optimization 
and controller design problems should be formulated and provides a good ground for 
further work. The completion of the recommended future work in the next section will 


make this approach to integrated plat/controller design a useful tool for engineers. 


29 




















B. RECOMMENDATIONS 


The convergence properties of the optimization algorithm are unclear and should 
be investigated. The first issue to consider is whether the algorithm converges to the 
global minimum for a convex problem. The second issue is to investigate the possi- 
bilities of finding a convex hull to the feasible set of plant parameters and controllers, 
and establish a lower bound for the optimal solution. 

The problem treated in this thesis only concerns the dynamics of the aircraft 
for small perturbations around a trimmed equilibrium. The formulation does not 
include the fact that the aircraft configuration given by the optimal value of ¢ must 
be possible to trim at the given flight condition. The problem formulation must be 
changed to include the requirement that the aircrft must be possible to trim. 

The performance requirements concerning maneuverability are translated to an 
H,. formulation based on the linearized model developed for small perturbations 
around an equilibrium. A non-linear dynamical model should instead be used for the 
maneurability requirements. 

The result of this work is a set of time invariant feedback controllers, each 
of them valid for one flight condition. This set of controllers can be used in the 
construction of a gain scheduled controller, that works for the non-linear aircraft 
dynamics at all flight conditions. The resulting closed loop system is time variant 
and non-linear. The gain scheduling must be done such that when this system is 
linearized at different flight conditions, the linear system has the desired properties. 
The linearized system dynamics should be as independent of the gain scheduling 
parameter as possible. A continuation of the work presented in this thesis is to find 


a way to interpolate between the different controller such that this is achieved. 
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APPENDIX A: MATLAB CODE 


In the appendix the Matlab code using LMI-Lab is listed. The main program is 
f4opt.m, which uses the two functions sol4.m and sol5.m. This combination of Matlab 
and LMI-Lab needs to be run on a SUN Sparc station and does not work on other 
computers. 


In f4opt.m the LMI’s are formulated using two equivalent versions of the LMI 





PAYEE YIAT(C) 4 BUQWEWTBIC) BIC) (CV! + DWH 
RW'*,Y',¢) = Bi (¢) =f 0 <0 
CY'+ DW? 0 —I 
these are 
_ AUQY'+ Y'AT(C) + B(QWi+ WTB (Ct 
RW, Y',¢) = +B; (¢)Bi7(¢) (CY'+ DW)? | <0 
CY’ + DW? —I 





that is used for solving for the controller, i.e. W and Y, and 


ae Ai(C)¥* + Y'A'T(¢) + By(Q)Wt + WT BSTC) + 
RW? Y',C) = +(C'y' + D'W)T (CY? + Diw') Bi(c) | <0 
Bi(¢)* —I 








that is used to solve for the optimal ¢. 
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4 
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This program performs an integrated minimization of elevator and spoiler, 
and calculation of H-infinity controller. The aircraft is an F-4 with 
aircraft data from Schmidt. The calculations are performed for five 
different flight conditions. 

The program uses the toolbox to be LMI-Lab and the two matlab functions 
sol4 and sol5. 

The system is modeled as 

x_dot=A*x+B_1*w+B_2*u 

z=C*x+D_2*u 

The program calculates the H-infinity controller K=Weinv(Y) (u=K*x) 

To use this program to design a controller the weights in the B_i, B_2, C and 
D_2 matrices can be changed. 

The output of the program is stored in the arrays lam {maximum eigen value 
for each LMI), zet (plant parameters) and cost (cost function). 


clear; 


g=32.174; % ft/s°2 
m=1210; % slugs 
c=16.0; % ft 


S= 


530; 4% ft~2 


Xdlc=-0.02*g*57.1/40; % drag from spoiler 
Zd1c=0.2*g*57.1/40; 4 destroyed lift from spoiler 
c=[1;1]; % weights in cost function 


rf 


ekkexx Flight condition 6 ***k*+ 


Iy=122190; % slug+ft 

U=1452; % ft/s 

Xu=-0.0072; % dimensional aerodynamic derivatives 
Xa=-29 .98; 

Zu=-0.515; 

Za=-716.7; 

Zad=-0.519; 

Zq=-2.24; 

Mu=0.0025; 

Ma=-28 .94; 

Mad=-0.122; 

Mq=-0.488; 

Xd=0.00; 

Zd=-70.67; 

Md=-16.00; 

xc=Md/Zd*Iy/m; % lever for elevator 
alpha0=2.6/57.3; % angle of attack in radians 
theta0=alpha0; 
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thrust_acc=0.05*g; % max acceleration from thrust 


sigma_g=1/57.3; % max commanded flight path angle 
sigma_a=5/U; % wind gusts 

umax1=20/57.3; 4% max elevator deflection 
umax2=40/57.3; 4% max spoiler deflection 
sigma_aout=1.5/57.3; % max error in alpha 


C=(zeros(3,6);0 1/sigma_aout 0 0 0 0;0 00 0 60 0;0 000 0 10]; 
D2=(1/umax1 0 0;0 i/umax2 0;0 O 1/thrust_acc;zeros(3,3)]; 


% A matrix for nominal plant 

An=[U*Xu Xa 0 -g*cos(theta0); U*Zu Za U+Zq -g*sin(theta0) ; 
UxMu Ma Mq 0; 0 0 1 0]; 

In=[U 0 0 0;0 U-Zad 0 0;0 -Mad 1 0; 0 0 0 1]; 

A006=([inv(In)*An zeros(4,2);0 10000;0 -10100]; 


% Calculation of derivatives for wing/body 
Za=Za-Zd; 

Ma=Ma-xc*Zdm/ly; 

Zq=Zq-xc*Zd/U; 

Mq=Mq-xc*xc*Zd¥m/Ty/U; 

Zad=Zad-xc*Zd/U; 

Mad=Mad-xc*xc*Zd*m/Iy/U; 


An=[U*Xu Xa 0 -g*cos(theta0); U*Zu Za U+Zq -g*sin(theta0) ; 
U*Mu Ma Mq 0; 0 0 1 0]; 

invIn=[1/U 0 0 0;0 1/(U-Zad) 0 0;0 Mad/(U-Zad) 10; 000 1]; 

A06=[invIn*An zeros(4,2);0 1000 0;0 -1010 0]; 


0 0; O Mad+xc*xc*m/Iy*(U-Zad) 0 0; 


1 0 0] /U/(U-Zad)~2; 
1 1/U 0; 0 xc#m/Iy xc¥*xc¥m/Iy/U 0;0 0 


0;0 0 
0;0 0]; 


0 
; 0 


? 


A16=[(dInd*Antinv(In)*dAnd)*Zd zeros(4,2); zeros(2,6)]; 


Bi06=zeros(6,2); Bi16=Bi06; B126=B106; 
B106(1:4,1)=A06(1:4,2)*sigma_a; 

B106(6,2)=sigma_g; 

aB106d=[0 0;1 0;xc¥*m/Iy 0;0 0]*sigma_a; 
B116=[(dInd*B1i06(1:4, :)+invIn*dB106d) *Zd;zeros(2,2)]; 


B206=[0 0 1/U;zeros(5,3)]; 
dB206d=[0 0 0;1 0 0;xc¥#m/Iy 0 0;0 0 0]; 
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dB206s=[0 Xdlc/Zdlc 0;0 1 0;0 0 0;0 0 O]; 
B216=[(dInd*B206(1:4, :)+invIn*dB206d)*Zd; zeros(2,3)]; 
B226=[invIn*dB206s*Zdlc; zeros(2,3)]; 

Jpeeeoeeoeorc Flight condition § seeeodooeoiiciiaaickicikic 


U=1167; 
Xu=-0.0136; 
Xa=-16.53; 
Zu=-0.747; 
Za=-848 .3; 
Zad=-1.24; 
Zq=-4.09; 
Mu=0.0022; 
Ma=-29.03; 
Mad=-0.288; 
Mq=-0.746; 
Xd=0.00; 
Zd=-90 .44; 
Md=-20.70; 
xc=Md/Zd*Iy/m; 
alphaO=1.6/57.3; 
theta0=alpha0; 
sigma_a=5/U; 


An=[U*Xu Xa 0 -g*cos(theta0); U*Zu Za U+Zq -g*sin(theta0) ; 
U*Mu Ma Mq 0; 0 0 1 0]; 

invIn=[1/U 0 0 0;0 1/(U-Zad) 0 0;0 Mad/(U-Zad) 1 0; 00 0 4]; 

AOS=[invIn*An zeros(4,2);0 10000;0 -10100]; 


dInd=[0 0 0 0;0 1 0 0; O Mad+xc*xc#m/Iy*(U-Zad) 0 0;0 0 0 0)/U/(U-Zad)~2; 
dAnd=[0 0 0 0;0 1 1/U 0; O xca#m/Iy xc*xc*m/Iy/U 0;0 0 0 0); 


Ai5=[(dInd*Ant+tinvIn*dAnd)*Zd zeros(4,2); zeros(2,6)]; 


B105=zeros(6,2); BiiS=Bi05; B1i25=Bi05; 
B105(1:4,1)=A05(1:4,2)*sigma_a; 

B105(6,2)=sigma_g; 

dBi05d=[0 0;1 0;xc#m/Iy 0;0 0]*sigma_a; 
B115=[(dInd*B105(1:4, :)+invIn*dB105d) *Zd;zeros(2,2)]; 


B205=[0 0 1/U;zeros(5,3)]; 

dB205d=[0 0 0;1 O 0;xc%m/Iy 0 0;0 0 0]; 

dB205s=[0 Xdlc/Zdlc 0;0 10;0 0 0;00 0]; 
B215=((dInd*B205 (1:4, :)+invIn*dB205d)*Zd; zeros(2,3)]; 
B225=([invIn*dB205s*Zdlc; zeros(2,3)]; 
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Joep oo ooOoOCee Flight condition 4 eee eol daa AAI IOI 


U=867 ; 
Xu=-0.0123; 
Xa=-7.43; 
Zu=-0.571; 
Za=-475 .4; 
Zad=-1.01; 
Zq=-2.89; 
Mu=-0 .0028; 
Ma=-7 .877; 
Mad=-0.234; 
Mq=-0 . 487; 
Xd=0.00; 
Zd=~49.65; 
Md=-11.40; 
xc=Md/Zd*Iy/m; 
alpha0=2.6/57.3; 
theta0=alpha0; 
sigma_a=5/U; 


An=[U*Xu Xa 0 -g*cos(theta0); U*Zu Za U+Zq -g*sin(theta0) ; 
UxMu Ma Mq 0; 0 0 1 0]; 
invIn=[1/U 0 0 0;0 1/(U-Zad) 
A04=[invIn*An zeros(4,2);0 1 


, 


0;0 Mad/(U-Zad) 10; 000 4); 
000;0 -10100)]; 


oo 


0 0] /U/(U-Zad)*2; 


dInd=[0 0 0 0; 0 
0 01; 


0 100; 0 Madtxc*xc#m/Iy*(U-Zad) 0 0;0 
dAnd=[0 0 0 0;0 1 1 0 


1 
1 1/U 0; 0 xc#m/Iy xc*xc¥m/Iy/U 0;0 
A14=[(dInd*An+invIn*dAnd)*Zd zeros(4,2); zeros(2,6)]; 


B104=zeros(6,2); Bi114=Bi04; B124=B104; 
B104(1:4,1)=A04(1:4,2)*sigma_a; 

B104(6,2)=sigma_g; 

dB104d=[0 0;1 0;xc*m/Iy 0;0 0]*sigma_a; 
B1i14=[(dInd+*B104(1:4, :)+invIn*dB104d)*Zd;zeros(2,2)]; 


B204=[0 0 1/U;zeros(5,3)]; 

dB204d=[0 0 0;1 0 0;xc*#m/Iy 0 0;0 0 0]; 

dB204s=[0 Xdlc/Zdlc 0;0 1 0;0 0 0;0 0 0]; 
B214=[(dInd*B204(1:4, :)+invIn*dB204d)*Zd; zeros(2,3)]; 
B224=([invIn#dB204s*Zdlc; zeros(2,3)]; 


JoesooROeeEeC Flight condition 3 eesgecidccogodiociccobiciiogctak 
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U=584; 
Xu=-0.0176; 
Xa=-24.25; 
Zu=-0.337; 
Za=-162.2; 
Zad=-0.591; 
Zq=-1.82; 
Mu=-0.0001; 
Ma=-1.927; 
Mad=-0.141; 
Mq=-0.307; 
Xd=0.00; 
Zd=-20 .98; 
Md=-4.90; 
xc=Md/Zd+Iy/m; 
alpha0=9.4/57.3; 
theta0=alpha0; 
sigma_a=5/U; 


An=[U*Xu Xa 0 -g*cos(theta0); U*Zu Za U+Zq -g*sin(theta0) ; 
U*Mu Ma Mq 0; 0 0 1 0]; 

invIn=[1/U 0 0 0;0 1/(U-Zad) 0 0;0 Mad/(U-Zad) 10; 000 1]; 

AO3=[invIn*An zeros(4,2);0 10000;0 -1010 0]; 


0 0; O Madtxc*xc*m/Iy*(U-Zad) 0 0;0 0 0 0]/U/(U-Zad)~2; 
1/U 0; 0 xc#m/Iy xc*xc*m/Iy/U 0;0 0 0 0]; 


. 


dInd=[0 0 0 0;0 1 
dAnd=[0 0 0 0;0 1 


A13=[(dInd*An+invIn*xdAnd)*Zd zeros(4,2); zeros(2,6)]; 


B1i03=zeros(6,2); B113=B103; B123=B103; 
Bi03(1:4,1)=A04(1:4,2)*sigma_a; 

B103(6,2)=sigma_g; 

dB103d=[0 0;1 0;xc*m/Iy 0;0 0]*sigma_a; 
B113=[(dInd*B103(1:4, :)+invIn*dB103d) *Zd;zeros(2,2)]; 


B203=[0 0 1/U;zeros(5,3)]; 

dB203d=[0 0 0;1 0 0;xc#m/Iy 0 0;0 0 0]; 

dB203s=[0 Xdlc/Zdlc 0;0 1 0;0 0 0;0 0 Oj; 
B213=[(dInd*B203(1:4, :)+invIn*dB203d)*Zd; zeros(2,3)]; 
B223=[invIn*dB203s*Zdlc; zeros(2,3)]; 


(aero oaodoaooodeo: Flight condition 2 saededkeioiko iar I 
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U=952; 
Xu=-0.0215; 
Xa=-§ .93; 
Zu=-1.170; 
Za=-1103.2; 
Zad=-2.00; 
Zq=-6 .00; 
Mu=-0 .0044; 
Ma=-17.00; 
Mad=-0.457; 
Mq=-0.993; 
Xd=0.00; 
Zd=-107.0; 
Md=-25.00; 
xc=Md/Zd*Iy/m; 
alpha0=0.5/57.3; 
theta0=alpha0; 
sigma_a=5/U; 


An=[U*Xu Xa 0 -g*cos(theta0); U*Zu Za U+Zq -g*sin(theta0) ; 
U*Mu Ma Mq 0; 0 0 1 0]; 

invIn=[1/U 0 0 0;0 1/(U-Zad) 0 0;0 Mad/(U-Zad) 10; 000 1]; 

AO2=[invIn*An zeros(4,2);0 10000;0 -10100]; 


dInd=[0 
dAnd=[0 


? 


0 00;0 1 0 0; O Mad+xc*xc*m/Iy*(U-Zad) 0 0;0 0 0 0]/U/(U-Zad)~2; 
000;0 1 1/U 0; 0 xc#m/Iy xc*xc*m/Iy/U 0;0 0 0 0]; 


A12=[(dInd*Ant+inviIn*dAnd)*Zd zeros(4,2); zeros(2,6)]; 


B102=zeros(6,2); B1i12=B102; B122=B102; 
B102(1:4,1)=A02(1:4,2)*sigma_a; 

B102(6,2)=sigma_g; 

dB1i02d=[0 0;1 0;xc#m/Iy 0;0 0]*sigma_a; 
Bi12=[(dInd*B102(1:4, :)+invIn*dB102d)*Zd;zeros(2,2)]; 


B202=[0 0 1/U;zeros(5,3)]; 

dB202d=[0 0 0;1 0 0;xc*m/Iy 0 0;0 0 0]; 

dB202s=[0 Xdlc/Zdlc 0;0 1 0;0 0 0;0 0 OJ; 
B212=[(dInd*B202(1:4,:)+invIn*dB202d)*Zd; zeros(2,3)]; 
B222=(invIn*dB202s*Zdlc; zeros(2,3)]; 

TTT T TIT ttt TTT Tt Tt TTT Lire ie te rer r Lelie rior rer orrr ere roc r er ses 
4 Calculations 

YEE OOO OO COG GOI OIG IO IOI IO OGIO GIO CIO IOI III IOI IIE 


¥% Inital value of zeta 











cost=[]; 
zet=(]; 


zeta=[0.1;0.7]; 
costold=c*zeta; 


for iter=1:40 


A2=A02+A12*zeta(1); 
B12=B102+B112*zeta(1); 
B22=B202+B212*zeta(1)+B222*zeta(2) ; 


A3=A03+A13*zeta(1); 
B13=B103+B113*zeta(1); 
B23=B203+B213*zeta(1)+B223*zeta(2) ; 


A4=A04+A14*zeta(1); 
B14=Bi04+B114*zeta(1); 
B24=B204+B214*zeta(1)+B224*zeta(2); 


AS5=A05+A15*zeta(1); 
B1i5=B105+B115*zeta(1); 
B25=B205+B215*zeta(1)+B225*zeta(2); 


A6=A06+A16*zeta(1); 
B16=B106+B116*zeta(1) ; 
B26=B206+B216*zeta(1)+B226*zeta(2) ; 


% Solve LMI’s for controller in iteration 1 sol 4 is used, whil sol5 is 
% used in the following iterations to correct the solution when feasp 

% gives a non-negative definite solution. In that case the old solution 
% is retained. 


if iter<2 

[Y2 W2 maxe2 x2]=sol4(A2,B12,B22,C,D2); 
[Y3 W3 maxe3 x3]=s014(A3,B13,B23,C,D2) ; 
[Y4 W4 maxe4 x4]=s014(A4,B14,B24,C,D2); 
{YS WS maxeS x5]=s014(A5,B15,B25,C,D2); 
CY6 W6 maxe6 x6]=s0l14(A6,B16,B26,C,D2) ; 
else 

(Y2 W2 maxe2 x2]=s015(A2,B12,B22,C,D2,x2); 
{Y3 W3 maxe3 x3]=s015(A3,B13,B23,C,D2,x3); 
[Y4 W4 maxe4 x4]=s015(A4,B14,B24,C,D2,x4); 
[YS WS maxeS x5J=s015(A5,B15,B25,C,D2,x5); 
CY6 W6 maxe6 x6]=s015(A6,B16,B26,C,D2,x6) ; 


38 








end; 


me=[maxe2 maxe3 maxe4 maxe5 maxe6]; 
lam=[lam; me]; 


%Solve LMI’s for optimal zeta 

R112=AO2*Y2+Y2*A02’ +B202*W2+W2’ *B202’ + (C*Y2+D2*W2) ? *(C*Y2+D2*W2) ; 
R113=A03*Y3+Y3%*A03’ +B2034*W3+W3’ *B203’ + (C¥Y3+D24*W3) ? *(C*Y3+D2*W3) ; 
R114=A04*Y4+Y4%A04’? +B2044*W4+W4’ *B204’ + (C¥Y¥4+D24W4) ? *(C*Y4+D24W4) ; 
R115=A05*Y5+Y5*A05’ +B205*WS+W5’ *B205’ + (C¥Y5+D24W5) ? *(C*Y5+D24W5S) ; 
R116=A06*Y6+Y6*A06’ +B206*W6+W6’ *B206’ + (C#Y6+D24W6) ? *(C#Y64+D2*W6) ; 


1m= [1 ; 

lm=add_lmi(1m,[6 2],8); 
lIm=add_lmi(1m, [6 2],8); 
lm=add_lmi(1m,[6 2],8); 
Im=add_lmi(im,(6 2] ,8); 
lm=add_lmi(l1m, [6 2],8); 


jim=add_lmi(1m, [1] ,1); 
lm=add_lmi(1m, [1] ,1); 


lm=add_var(1m, 1, [1 1]); 
lm=add_var(im, 1, [1 1]); 


im=add_term(1m, [1 1 1 0], R112); 
Im=add_term(1m, [1 1 1 1], B212, W2, ’s’); 
lm=add_term(1m, [1 11 1], A12, Y2, ’s’); 
lm=add_term(1m, [1 1 1 2], B222, W2, ’s’); 
lm=add_term(im, [1 1 2 0], B102); 
lm=add_term(1m, [1 1 2 1], B112); 
Im=add_term(1m, [1 2 2 0], -1); 
lm=add_term(1m, [2 1 1 0], R113); 
jJm=add_term(1m, [2 1 1 1], B213, W3, ’s’); 
lm=add_term(1m, [2 11 1], A13, Y3, ’s’); 
Im=add_term(1m, [2 1 1 2], B223, W3, ’s’); 
Im=add_term(1m, [2 1 2 0], B103); 
lm=add_term(1m, [2 1 2 1], B113); 
Im=add_term(1m, [2 2 2 0], -1); 
lm=add_term(1Jm, [3 1 1 0], R114); 
Im=add_term(1m, (3 1 1 1], B214, W4, ’s’); 
lm=add_term(Im, [3 1 1 1], A14, Y4, ’s’); 
Im=add_term(1m, [3 1 1 2], B224, W4, ’s’); 
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lm=add_term(1m, 
lim=add_term(1m, 
lm=add_term(1m, 


lm=add_term(1m, 
lm=add_term(1m, 
lIm=add_term(1m, 
lm=add_term(1m, 
lm=add_term(ln, 
Im=add_term(lm, 
lm=add_term(1m, 


lm=add_term(1m, 
lm=add_term(1m, 
im=add_term(l1m, 
lm=add_term(lm, 
lm=add_term(1m, 
lm=add_term(1m, 
lm=add_term(1m, 


lm=add_term(1m, 
lm=add_term(im, 


NO ooo en oe on 


NON YN FP FP eS eS 


NO NM NY BY RR 





B104); 
Bi14); 
-1); 


R115); 

B215, W5, ’s’); 
A1S, YS, ’s’); 

B225, WS, ’s’); 
B1i05); 

B1ii5); 

-1); 


R116); 

B2i6, W6, ’s’); 
A1i6, Y6, ’s’); 
B226, W6, ’s’); 
B106) ; 

B1i16); 

-1); 


[copt ,xopt]=linobj(1m,c,[0 0 0 50 1],zeta); 


iter 
zeta 


%Controllers for each flight condition 


K2n=W2*inv(Y2); 
K3n=W3*inv(Y3) ; 
K4n=W4*inv(Y4) ; 
K5n=WS*inv(Y5) ; 
K6n=W6*inv(Y6) ; 
K2=[K2; K2n]; 
K3=[K3; K3n]; 
K4=(K4; K4n]; 
K5=[K5; K5n]; 
K6=[K6; K6n]; 


Cost=c*xopt 
cost=[cost; Cost]; 
zet=[zet; zeta’]; 
xopt 
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zeta=xopt ; 


if abs (Cost-Costold) /Costold<0.0001 
if max(me)<0 

break; 

end; 

end; 


% Closed loop system in the final version 
Ac2=A2+B22*K2n ; 
Ac3=A3+B23*K3n ; 
Ac4=A4+B24*K4n ; 
AcS=A5+B25*K5n; 
Ac6=A6+B26*K6n; 


damp (Ac2) 
damp (Ac3) 
damp (Ac4) 
damp (Ac5) 
damp (Ac6) 
D=zeros (3,3); 


end 


4] 











function [Y, W, maxe, x]=sol14(A,B1,B2,C,D) 
% Solves the H-inf control problem and returns Y, W matrices 
% the maximum eigenvalue and the solution vector in LMI format 
im=[] ; 

lm=add_lmi(1m, [6 6] ,12); 

1m=add_1mi(1m, [6] ,6) ; 

lm=add_var(1m, 1, [6 1]); 

Im=add_var(1lm, 2, [3 6]); 

lm=add_term(1m,[1 1 1 0], B1*B1’); 
lm=add_term(1m, [1 1 1 1], A, eye(6), ’s’); 
lm=add_term(1m, [1 1 1 2], B2, eye(6), ’s’); 
lm=add_term(1im, [1 2 1 1], C); 
lm=add_term(1m, [1 2 1 2], D); 
lm=add_term(1m, [1 2 2 0], -1); 
lm=add_term(lm, [-2 1 1 1], eye(6), eye(6)); 
{tmin,xfeas]=feasp(1m,[200 0 100 1],0); 
evalsys=eval_1mi(1n, [xfeas]) ; 
(lhs,rhs]=show_lmi(evalsys, 1); 

Y=dec2mat (1m,xfeas,1); 

W=dec2mat (1m,xfeas ,2) ; 

maxe=max (eig(lhs)) 

x=xfeas ; 
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function [Y, W, maxe, x]=so015(A,B1,B2,C,D,x1) 

% Solves the H-inf control problem given A, Bi, B2, C, D matrices 

% and the previous LMI solution. This function checks if all the 

% eigenvalues are negative and if they are not returns a solution based 


% on the previous iteration. 
1m=(] ; 

lm=add_lmi(1m, [6 6] ,12); 
lm=add_1mi(1m, [6] ,6); 
lm=add_var(lm, 1, [6 1]); 
lm=add_var(1m, 2, [3 6]); 
lm=add_term(1m,[1 11 0], B1*B1’); 


lm=add_term(lm, [1 11 1], A, eye(6), ’s’); 
lm=add_term(1m, [1 1 1 2], B2, eye(6), ’s’); 


lm=add_term(lm, [1 2 1 1], ¢); 
lm=add_term(1m, [1 2 1 2], D); 
lm=add_term(1m, [1 2 2 0], -1); 


lm=add_term(lm, [-2 1 1 1], eye(6), eye(6)); 


evalsys=eval_l1mi (1m, [x1]); 
(lhs,rhs]=show_lmi(evalsys, 1); 
maxe1=max(eig(lhs) ) 


([tmin,xfeas]=feasp(1m, [200 0 100 1],0); 


evalsys=eval_1mi(1m, [xfeas]) ; 
[lhs ,rhs]=show_lmi(evalsys, 1); 
maxe=max (eig(lhs) ) 


if maxe<0O 
x=xfeas; 
else 
maxe=maxel ; 
x=x1i; 

end; 


Y=dec2mat(1m,x,1); 
W=dec2mat (1m,x,2); 
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